home *** CD-ROM | disk | FTP | other *** search
/ Mac100% 1998 November / MAC100-1998-11.ISO.7z / MAC100-1998-11.ISO / ｵﾝﾗｲﾝｿﾌﾄ定点観測 / ﾕｰﾃｨﾘﾃｨ / Mops 3.2.sea / Mops 3.2 / Mops source / PPC source / zFloating point < prev    next >
Text File  |  1998-04-29  |  19KB  |  823 lines
  1. ¥  Mar 98  mrh  Initial PowerPC version.
  2.  
  3. :f fNum?  false  ;f
  4.  
  5. syscall    dec2num
  6. syscall    num2dec
  7. syscall    dec2str
  8. syscall    str2dec
  9.  
  10. $ BD42        simple_op    SF@
  11. $ BD43        simple_op    SF!
  12.  
  13. 8  constant    1FLOAT
  14.  
  15. ¥ Our default float is an IEEE double.  Thus DF@ etc. are aliases for F@ etc.
  16.  
  17. : DF@        inline{ f@}  ;
  18. : DF!        inline{ f!}  ;
  19.  
  20. ¥ SF@ and SF! are in pnuc1.
  21.  
  22. : FLOAT+    inline{ 8 +} ;
  23. : FLOAT-    inline{ 8 -} ;
  24. : FLOATS    inline{ 3 <<} ;
  25. : FALIGN    align8  ;
  26. : FALIGNED    #align8  ;
  27.  
  28. : DFLOAT+    inline{ 8 +} ;
  29. : DFLOAT-    inline{ 8 -} ;
  30. : DFLOATS    inline{ 3 <<} ;
  31. : DFALIGN    align8  ;
  32. : DFALIGNED    #align8  ;
  33.  
  34. : SFLOAT+    inline{ 4 +} ;
  35. : SFLOAT-    inline{ 4 -} ;
  36. : SFLOATS    inline{ 2 <<} ;
  37. : SFALIGN    align4  ;
  38. : SFALIGNED    #align4  ;
  39.  
  40.  
  41. : F,  ( F: r -- )
  42.     falign  DP f!  1float ++> DP  ;
  43.  
  44. : FVARIABLE
  45.     falign  variable  ;
  46.  
  47.  
  48. ¥            =====================================
  49. ¥                           FP mode setting
  50. ¥            =====================================
  51.  
  52. ¥ rounding modes:
  53. enum{  roundToNearest  roundToZero  roundUp  roundDown  }
  54.  
  55. :ppc_code (fpscr)
  56.     fr0                mffs,
  57.     fr0        0  r4    stfd,
  58.                     blr,
  59. ;ppc_code
  60.  
  61. :ppc_code (>fpscr)
  62.     fr0        0    r4    lfd,
  63.     $ FF     fr0        mtfsf,
  64.                     blr,
  65. ;ppc_code
  66.  
  67. : FPSCR
  68.     ftemp  (fpscr)  4+ @
  69. ;
  70.  
  71. : >FPSCR
  72.     ftemp 4+ !  ftemp (>fpscr)  drop  ;
  73.  
  74. : rounding_mode
  75.     fpscr  3 and  ;
  76.  
  77. : >rounding_mode
  78.     3 and fpscr $ ffffffc0 and or  >fpscr  ;
  79.  
  80.  
  81.  
  82. ¥            =====================================
  83. ¥                        Rounding to integer
  84. ¥            =====================================
  85.  
  86. variable    2**52    4 reserve
  87. variable    -2**52    4 reserve
  88.  
  89. $ 43300000  2**52  !
  90. $ C3300000    -2**52 !
  91.  
  92. : round  { %x ¥ %bias -- %x' }
  93.     %x f0>=
  94.     IF                ¥ positive - add 2**52, then subtract it
  95.         2**52 f@ -> %bias
  96.         %x %bias f>
  97.         IF    %x  EXIT  THEN
  98.     ELSE
  99.         -2**52 f@  -> %bias
  100.         %x %bias f<
  101.         IF  %x  EXIT  THEN
  102.     THEN
  103.     %x %bias f+  %bias f-
  104. ;
  105.  
  106. : FLOOR
  107.     fpscr
  108.     roundDown >rounding_mode
  109.     round
  110.     >fpscr
  111. ;
  112.  
  113. : FROUND
  114.     fpscr
  115.     roundToNearest >rounding_mode
  116.     round
  117.     >fpscr
  118. ;
  119.  
  120.  
  121. ¥            =====================================
  122. ¥                   FP to/from decimal conversion
  123. ¥            =====================================
  124.  
  125. ¥ Some useful constants:
  126.  
  127.  256    constant    NEG
  128.    0    constant    POS
  129.  
  130. 1    constant    fixedDecimal
  131. 0    constant    floatDecimal
  132.  
  133.  
  134. ¥ Our DEC class combines the SANE/MathLib decimal and decform records,
  135. ¥  and adds appropriate methods to do what we want.  We could equally
  136. ¥  well have made these two separate classes, but this is a little
  137. ¥  less complicated, especially as we only instantiate one Dec object.
  138.  
  139. :class  DEC super{ object }
  140. public
  141.     ¥ decimal record  ( x = (-1)^sgn * 10^exp * digits )
  142. record
  143. {
  144.     byte    sign
  145.     byte    unused
  146.     int        exp            ¥ offs 2
  147.     ubyte    length        ¥ offs 4
  148. 36    bytes    text        ¥ PPC Numerics manual sez 36 is the max len
  149.     ubyte    unused1
  150.     uint    index        ¥ offs 42 - used by the scanner
  151.     
  152.     ¥ decform record starts here (we're aligned, at offs 44)
  153.  
  154.     byte    style
  155.     byte    unused2
  156.     int        #digits        ¥ # of sig digits, if float;
  157.                         ¥ # dec. places, if fixed.
  158.     int        valid?        ¥ used by the scanner.
  159. }
  160.  
  161.  
  162. :m CLEAR:
  163.     addr: sign 41 erase  ;m
  164.  
  165. :m EINIT:    clear: self  FloatDecimal  put: style  19  put: #digits  ;m
  166. :m FINIT:    clear: self  FixedDecimal  put: style  ;m
  167.  
  168. :m SETSTYLE:    put: style    ;m
  169. :m SET#DIGITS:    put: #digits  ;m
  170. :m SETEXP:        put: exp    ;m
  171. :m EXP:            get: exp    ;m
  172. :m SIGN:        get: sign    ;m
  173.  
  174. :m ZERO:        ¥ Puts a "0" in decimal record (length 1)
  175.     clear: self  $ 0130  addr: length  w!  ;m
  176.  
  177. :m FLOAT:  ( -- d )
  178.     ^base                        ¥ Addr of decimal record
  179.     dec2num                        ¥ returns a double
  180. ;m
  181.  
  182. ¥ >DEC: converts the passed-in double float to decimal.
  183.  
  184. :m >DEC:
  185.     addr: style                ¥ addr of decform record
  186.                             ¥ (double is on top of FP stack)
  187.     ^base                    ¥ Addr of decimal record
  188.     num2dec
  189. ;m
  190.  
  191. ¥ Ascii input
  192.  
  193. :m $>DEC:        ¥ ( addr len -- )
  194.     str255  1+
  195.     clear: index  addr: index
  196.     ^base  addr: valid?
  197.     str2dec  ;m
  198.  
  199.  
  200. :m $>DEC?:  { addr len -- b }
  201.         ¥ Attempts to convert the passed-in string, using $>DEC:.
  202.         ¥ Returns True if all the input was read.  Otherwise
  203.         ¥ we assume the terminating (non-scanned) character is
  204.         ¥ invalid, and return False.
  205.     addr len  $>dec: self
  206.     get: valid? 0<>            ¥ turn a C flag into a Mops one
  207. ;m
  208.  
  209. ¥ Ascii output
  210.  
  211. :m FORMAT:    { ¥ addr ^c -- addr' len }
  212.     pad -> addr        ¥ as good a place as any for the output
  213.     addr: style        ¥ addr of decform
  214.     ^base            ¥ addr of decimal record
  215.     addr
  216.     dec2str
  217.     addr -> ^c
  218.     BEGIN  ^c c@  WHILE  1 ++> ^c  REPEAT
  219.     addr  ^c addr -  ;m
  220.  
  221. :m PRINT:
  222.     format: self  type  ;m
  223.  
  224. :m DUMP:
  225.     ." sign:    "    get: sign  IF  & -  ELSE  & +  THEN  emit  cr
  226.     ." exp:     "    print: exp        cr
  227.     addr: length  count  type  cr
  228.     ." style:   "    get: style  IF  ." fixed"  ELSE  ." float"  THEN  cr
  229.     ." index:   "    get: index  .  cr
  230.     ." #digits: "    get: #digits  .  cr  ;m
  231.     
  232. ;class
  233.  
  234. dec    theDec
  235.  
  236. : #DIGITS    set#digits: theDec  ;
  237.  
  238. : E.R  ( F: r -- ) { wid ¥ svOut -- }
  239.     out -> svOut
  240.     floatDecimal  setStyle: theDec
  241.     wid 6 -  #digits        ¥ Allow for point, sign, and e+nn
  242.     >dec: theDec
  243.     print: theDec
  244.     wid  out svOut -  -  spaces  ;
  245.  
  246. : E.  ( F: r -- )    23 e.r  ;        ¥ 64-bit double gives us a 52-bit
  247.                                     ¥  fraction, or 17 decimal digits of
  248.                                     ¥  significance.
  249.                                     
  250. : FS.    inline{ e.}  ;                ¥ ANSI synonym.
  251.  
  252.  
  253. : F.R  ( F: r -- )  { wid ¥ #dig svOut -- }
  254.     out -> svOut
  255.     floatDecimal  setStyle: theDec
  256.     wid 2-  #digits            ¥ Allow for sign and dec point
  257.     >dec: theDec
  258.     fixedDecimal setStyle: theDec
  259.     exp: theDec  negate  dup -> #dig  #digits
  260.     sign: theDec NIF  space  THEN
  261.     #dig NIF  space  THEN        ¥ In this case, no dec point
  262.     print: theDec
  263.     wid  out svOut -  -  spaces  ;
  264.  
  265.  
  266. : F.  ( F: r -- )    ¥ Default fixed-point format print.
  267.     19 f.r  ;        ¥ 17 digits of precision, plus the "." and the
  268.                     ¥  trailing space gives a field width of 19.
  269.  
  270. : REPRESENT  { addr len %x -- exp sign ok? }
  271.         ¥ ANSI word to give a generic ASCII representation for
  272.         ¥  a floating point number.  The fraction, rounded to len
  273.         ¥  digits, is placed at addr (using round to nearest rounding).
  274.         ¥ Returns the exponent (assuming the decimal point is just to
  275.         ¥  the left of the fraction), the sign (true = minus) and
  276.         ¥  a success/failure flag.
  277.     fpscr
  278.     roundToNearest >rounding_mode
  279.     len #digits
  280.     %x >dec: theDec
  281.     >fpscr
  282.     exp: theDec  get: ivar> length in theDec +
  283.     sign: theDec
  284.     addr: ivar> text in theDec  dup c@
  285.     & 0  & 9  within? nip
  286.     IF        addr len  cmove    true
  287.     ELSE    drop  addr len  & * fill    false
  288.     THEN
  289. ;
  290.  
  291.  
  292.  
  293. (*    >FLOAT Attempts to convert the passed-in ASCII string to
  294.     floating, if possible.  On success, returns true on the data
  295.     stack and the floating number on the FP stack.  On failure, 
  296.     returns false on the data stack.
  297.     
  298.     The ANSI standard defines this word to be very liberal - I think
  299.     it's more liberal than str2dec in MathLib.  So eventually we
  300.     might have to add a preprocessing scan on the input string.
  301.     But as it stands it will work for any reasonable-looking number.
  302. *)
  303.  
  304. : >FLOAT  { addr len -- b }  ( F: -- r  | -- )
  305.                 ¥ Converts the passed-in ASCII string to
  306.                 ¥ floating, if possible.  On success, returns
  307.                 ¥ true on the data stack and the floating number
  308.                 ¥ on the FP stack.  On failure, returns false
  309.                 ¥ on the data stack.
  310.  
  311.     addr len  $>dec?: thedec  NIF  false  EXIT  THEN
  312.     float: theDec  true  ;
  313.  
  314.  
  315. ¥            ==============================
  316. ¥                      Interpretation
  317. ¥            ==============================
  318.  
  319.  
  320. : FNUMBER    ¥ ( addr -- flt T  |  -- F )
  321.             ¥ Attempts to convert token at addr to a float.
  322.     count  >float  ;
  323.  
  324.  
  325.  
  326. : FLITERAL  ( F: r -- )            ¥ Compiles an in-line float.
  327.     ftemp f!                    ¥ store float at ftemp
  328.     const_data_ref  postpone f@    ¥ compile fetch from const_data area
  329.     ftemp 8  add: const_data    ¥ store the float in this defn's
  330.                                 ¥  const_data
  331. ;        immediate
  332.  
  333.  
  334. (*    FNUM? is called from INTERPRET, to check if the string at addr
  335.     is a FP number.  ANSI specifies that it's a float if it contains
  336.     "E" and the base is decimal.  As well as this, Mops has always 
  337.     taken an embedded decimal point to indicate a float, if floating 
  338.     point is loaded.  So we do all these here.  We assume it's a
  339.     float if we're in decimal, and the string contains E, e, or a
  340.     decimal point.
  341.     If these conditions are met but fNumber returns an error, we
  342.     give a "not found" error, since at this point INTERPRET has already
  343.     looked in the dictionary for the string and not found it, and
  344.     if it contains E or . then it certainly isn't a legal integer.
  345. *)
  346.  
  347. :f FNUM?  { addr ¥ addr' len -- T  |  -- addr F } ( F: -- r  | -- )
  348.  
  349.     base 10 <> IF  addr  false  EXIT  THEN
  350.  
  351.     addr count  -> len  -> addr'
  352.     addr' len  & .  scan  nip
  353.     NIF    addr' len  & E  scan  nip
  354.         NIF    addr' len  & e  scan  nip
  355.             NIF    
  356.                 addr  false  EXIT
  357.             THEN
  358.         THEN
  359.     THEN
  360.  
  361.     addr fnumber  ?notFound
  362.     state  IF  postpone fliteral  THEN  true  ;f
  363.  
  364.  
  365. ¥ FP stack dump:
  366.  
  367. : dummy ;
  368.  
  369. : F.S  { ¥ start-addr end-addr ¥ val dpth -- }
  370.  
  371.     0.0 0.0  dummy        ¥ push off the 2 cells in regs
  372.     0 -> out
  373.     fsp -> start-addr  fsp0 -> end-addr
  374.  
  375.     start-addr end-addr >
  376.     IF    ."  underflow"    cr f2drop  EXIT  THEN
  377.     start-addr end-addr =
  378.     IF    ."  empty"        cr f2drop  EXIT  THEN
  379.  
  380.     end-addr start-addr - 3 a>>  -> dpth
  381.     ."  depth "  dpth .
  382.  
  383.     BEGIN
  384. ¥        ?pause
  385.         cr
  386.         start-addr f@ e.
  387.  
  388.         8 ++> start-addr
  389.         start-addr end-addr >=
  390.     UNTIL
  391.     cr  f2drop
  392. ;
  393.  
  394. ¥ F.S+ does our FP stack display in the Mops window.  It's forward
  395. ¥ defined, and initially resolved to a dummy word which does nothing.
  396. ¥ Here we re-resolve it.  It's called from the DS: method in 
  397. ¥ class TEFwind.
  398.  
  399. :f F.S+
  400.     getbotx: tempRect 2/ negate  0  setOrigin
  401.     10 10 gotoxy  ." FP stack:  "
  402.     f.s
  403.     0 0 setOrigin
  404. ;f
  405.  
  406. ¥            ===========================================
  407. ¥                    Integer <-> FP conversions
  408. ¥            ===========================================
  409.  
  410. (*    The algorithm here looks a bit magic, but isn't really complicated.
  411.     The initial problem is the different way we represent negative
  412.     numbers in FP and integer arithmetic.  For a single (4-byte) integer,
  413.     we first bias the integer by flipping the top bit.  This converts 
  414.     the 2's complement integer to an an unsigned integer biassed by 
  415.     2**31.
  416.     
  417.     We now convert this to FP by prepending another 32-bit word
  418.     which contains just the exponent.  To get this exponent, remember
  419.     that if the binary point is straight after the exponent (just
  420.     before bit 12) the biassed exponent is 1023 (see the manual).
  421.     So moving the binary point down to the integer position, before
  422.     the (imaginary) bit 64, we get an exponent of
  423.     1023+64-12 = 1024+51 = $433.
  424.  
  425.     Now the new FP number has an implicit leftmost 1 in the (52-bit)
  426.     fraction part, so in doing this we've in effect added another
  427.     bias of 2**52.  All we then have to do is load this into an FPR,
  428.     and subtract the total bias, which at the same time normalizes 
  429.     the result for us.  The total bias is 2**52 + 2**31.  In FP,
  430.     this is a number with the same exponent as above ($433), and
  431.     a 1 in bit 32 (the 2**52 is again implicit).  This is the number
  432.     "fmagic" below.
  433.     
  434.     For converting a double integer to floating, we do the equivalent.
  435.     THe Forth Standard says it's ambiguous if the number can't be
  436.     exactly represented as a float, which for us means if there are
  437.     more than 52 significant bits in the double integer.  In this case,
  438.     as it's ambiguous, we can do anything we like.  Our algorithm will
  439.     simply give a rounded result, which is pretty reasonable.
  440. *)
  441.  
  442.         variable  fmagic    8 allot
  443. $ 43300000    fmagic !
  444. $ 80000000    fmagic 4+ !
  445.  
  446.         variable  fmagic2    8 allot
  447. $ 45300000    fmagic2 !
  448. $ 80000000    fmagic2 4+ !
  449.  
  450.         variable  fmagic3    8 allot
  451. $ 43300000    fmagic3 !
  452. 0            fmagic3 4+ !
  453.  
  454. : S>F  ( n -- )  ( f: -- r )
  455.     $ 43300000  ftemp !            ¥ Store exponent
  456.     $ 80000000  xor                ¥ flip top bit of integer
  457.     ftemp 4+ !                    ¥ store as lo half of FP number
  458.     ftemp f@                    ¥ fetch to FP regs
  459.     fmagic f@  f-                ¥ subtract FP bias
  460. ;
  461.  
  462. : D>F
  463.     $ 45300000  ftemp !            ¥ Store exponent
  464.     $ 80000000  xor                ¥ flip top bit of integer
  465.     ftemp 4+ !                    ¥ store as lo half of FP number
  466.     ftemp f@                    ¥ fetch to FP regs
  467.     fmagic2 f@  f-                ¥ subtract FP bias
  468.     
  469.     $ 43300000  ftemp !            ¥ Store exponent
  470.     ftemp 4+ !                    ¥ store as lo half of FP number
  471.     ftemp f@                    ¥ fetch to FP regs
  472.     fmagic3 f@  f-                ¥ subtract unsigned FP bias
  473.     
  474.     f+                            ¥ add both components
  475. ;
  476.  
  477.  
  478. ¥ F>D could be implemented along the lines of ROUND (see above), but 
  479. ¥  with all the monkeying around with rounding modes etc, it turns out
  480. ¥  to take 4 times as many instructions as the following, which we've
  481. ¥  lifted from Metrowerks:
  482.  
  483.  
  484. :ppc_code F>D
  485.  
  486.     fr2        -8        rRP        stfd,
  487.     fr2        fr1                fmr,
  488.     fr1        0        rFSP    lfd,
  489.     rFSP    rFSP    8        addi,
  490.  
  491.     r3        -4        rSP        stw,
  492.     r4        -8        rSP        stwu,
  493.  
  494.     r4        -8        rRP        lwz,
  495.     r3        -4        rRP        lwz,
  496.     r5        r4 12 21 31        rlwinm,
  497.     r5        1023            cmpli,
  498. lt if,
  499.     r3        0                li,
  500.     r4        0                li,
  501.                             blr,
  502. then,
  503.  
  504.     r6        r4                mr,
  505.     r4        r4 0 12 31        rlwinm,        ¥ mantissa has implied 1 bit
  506.     r4        r4        $ 10    oris,        ¥ [nh = (nh & 0xfffff) | 0x100000]
  507.  
  508. ¥ when exp-1075 < 0 we need to shift right, else shift left
  509.  
  510.     r5        r5        -1075    addi,
  511.     r5        0                cmpi,
  512. lt if,
  513. ¥ shift r3:r4 right by -exp bits
  514.  
  515.     r5        r5                neg,
  516.     r8        r5        32        subfic,
  517.     r9        r5        -32        addic,
  518.     r3        r3        r5        srw,        
  519.     r10        r4        r8        slw,
  520.     r3        r3        r10        or,    
  521.     r10        r4        r9        srw,    
  522.     r3        r3        r10        or,        ¥ low word
  523.     r4        r4        r5        srw,    ¥ high word
  524.  
  525. else,
  526. ¥ shift r3:r4 left by 'exp' bits
  527. ¥ if exp >= 32, then the result is undefined. later we might force
  528. ¥ this to the maximum long long value or whatever might make sense
  529. ¥ to users.
  530.     r8        r5        32        subfic,
  531.     r9        r5        -32        addic,
  532.     r4        r4        r5        slw,
  533.     r10        r3        r8        srw,
  534.     r4        r4        r10        or,
  535.     r10        r3        r9        slw,
  536.     r4        r4        r10        or,            ¥ high word
  537.     r3        r3        r5        slw,        ¥ low word
  538. then,
  539.  
  540. ¥ if the sign bit was set, negate the long long
  541.     r6        r6    0 0 0        rlwinm.,
  542. eq bclr,
  543.  
  544.     r3        r3        0        subfic,
  545.     r4        r4                subfze,
  546.  
  547.                             blr,
  548. ;ppc_code
  549.  
  550. : F>S    f>d  drop  ;
  551.  
  552.  
  553. ¥                ==================================
  554. ¥                              Constants
  555. ¥                ==================================
  556.  
  557. ¥ 0.0 is already defined - we keep it in the register fr14.
  558.  
  559. 1.0 fconstant    F1.0        ¥ FSL uses this
  560.  
  561.  
  562. ¥                ==================================
  563. ¥                          Transcendentals etc.
  564. ¥                ==================================
  565.  
  566. syscall sin
  567. syscall cos
  568. syscall tan
  569.  
  570. syscall    asin
  571. syscall    acos
  572. syscall    atan
  573. syscall atan2
  574.  
  575. syscall    sinh
  576. syscall    cosh
  577. syscall tanh
  578.  
  579. syscall    asinh
  580. syscall acosh
  581. syscall    atanh
  582.  
  583. syscall    exp
  584. syscall    expm1
  585. syscall    pow
  586. syscall    log
  587. syscall    log1p
  588. syscall    log10
  589.  
  590. syscall    sqrt
  591.  
  592.  
  593. : FSIN        sin  ;
  594. : FCOS        cos  ;
  595. : FTAN        tan  ;
  596.  
  597. : FASIN        asin  ;
  598. : FACOS        acos  ;
  599. : FATAN        atan  ;
  600. : FATAN2    atan2 ;
  601.  
  602. : FSINH        sinh  ;
  603. : FCOSH        cosh  ;
  604. : FTANH        tanh  ;
  605.  
  606. : FASINH    asinh  ;
  607. : FACOSH    acosh  ;
  608. : FATANH    atanh  ;
  609.  
  610. : FEXP        exp        ;
  611. : FEXPM1    expm1    ;
  612. : F**        pow        ;
  613. : FLN        log        ;        ¥ log to base e
  614. : FLNP1        log1p    ;
  615. : FLOG        log10    ;        ¥ log to base 10
  616.  
  617. : FSQRT        sqrt    ;
  618.  
  619. : 1/F        1.0 fswap f/  ;
  620. : F**2        inline{ fdup f*}  ;
  621.  
  622.  
  623.  
  624. ¥            ======================================
  625. ¥                    Infinities and NANs
  626. ¥            ======================================
  627.  
  628. 1.0 0.0 f/            fconstant    infinity
  629. infinity fnegate    fconstant    -infinity
  630.  
  631.  
  632. : SNAN    ( code -- ) ( F: -- NAN<code> )
  633.     8 <<  $ 7ff00000 or  ftemp !  ftemp sf@  ;
  634.  
  635. : QNAN    ( code -- ) ( F: -- NAN<code> )
  636.     8 <<  $ 7ff80000 or  ftemp !  ftemp sf@  ;
  637.  
  638. 100 sNAN    fconstant    UNDEF
  639.  
  640.  
  641. ¥            ===================================
  642. ¥                        Class Float
  643. ¥            ===================================
  644.  
  645. (*    Class Float allows a floating value to be a high-level object, which
  646.     means it can be an ivar or array element.  In 68k Mops we included
  647.     methods like *: and /:, but as these involve storing the updated
  648.     float to memory, which is much slower on the PowerPC than doing the
  649.     operations, we're now omitting them.  We just include +: and -:
  650.     for consistency with our other numeric types.
  651. *)
  652.  
  653.  
  654. :class  FLOAT  super{ object }
  655.  
  656. 8    bytes    data
  657.  
  658. :m GET:        ¥ ( -- x )    Pushes private data onto FP stack
  659.     inline{ ^base f@}  ;m
  660.  
  661. :m PUT:        ¥ ( x -- )  Stores float into private data
  662.     inline{ ^base f!}  ;m
  663.  
  664. :m ->:        ¥ ( float -- )  Assigns value of passed-in Float to this Float
  665.     inline{ f@ ^base f!}  ;m
  666.  
  667.  
  668. ¥ ----- Arithmetic operations take a stack float (not a Float obj)
  669.  
  670. :m +:
  671.     inline{ ^base f@ f+ ^base f!}  ;m
  672.  
  673. :m -:
  674.     inline{ ^base f@ f- ^base f!}  ;m
  675.  
  676. :m ABSVAL:    ¥ ( -- abs )  Returns absolute value.
  677.     inline{ ^base f@ fabs}  ;m
  678.  
  679. :m ABS:        ¥ ( -- )  Replaces obj's data with its absolute. Doesn't
  680.                 ¥            return anything.
  681.     inline{ ^base f@ fabs ^base f!}  ;m
  682.  
  683. :m NEG:        ¥ ( -- val )  Returns object value with sign negated
  684.     inline{ ^base f@ fnegate}  ;m
  685.  
  686. :m  NEGATE:    ¥ ( -- )  Negates the object's data. Doesn't return anything.
  687.     inline{ ^base f@ fnegate ^base f!}  ;m
  688.  
  689. :m PRINT:    ^base f@  e.  ;m
  690.  
  691. ;class
  692.  
  693.  
  694. ¥            =================================
  695. ¥                    Floating arrays
  696. ¥            =================================
  697.  
  698. :class     FARRAY    super{ indexed-obj }  1float indexed
  699.  
  700. :m AT:  ( index -- n )        inline{ ^elem f@}  ;m
  701. :m TO:  ( n index -- )        inline{ ^elem f!}  ;m
  702. :m +TO:  ( n index -- )        inline{ ^elem dup f@ f+ f!}  ;m
  703. :m -TO:  ( n index -- )        inline{ ^elem dup f@ f- f!}  ;m
  704.  
  705.  
  706. :m FILL:        ¥ ( value -- )  Fills all elements with value.
  707.     idxbase  limit floats  bounds
  708.     ?DO  fdup  i f!  1float +LOOP  fdrop  ;m
  709.  
  710. :m PRINT:    ¥ Prints all elements
  711.     limit: self 0 ?DO  i dup  4 .r  space  at: self  e. cr
  712.     LOOP  ;m
  713.  
  714. :m CLASSINIT:
  715.     undef
  716.     limit: self  FOR  fdup  i to: self  NEXT  fdrop  ;m
  717.  
  718. ;class
  719.  
  720.  
  721. :class     SFARRAY    super{ indexed-obj }  4 indexed
  722.  
  723. :m AT:  ( index -- n )        inline{ ^elem sf@}  ;m
  724. :m TO:  ( n index -- )        inline{ ^elem sf!}  ;m
  725. :m +TO:  ( n index -- )        inline{ ^elem dup sf@ f+ sf!}  ;m
  726. :m -TO:  ( n index -- )        inline{ ^elem dup sf@ f- sf!}  ;m
  727.  
  728.  
  729. :m FILL:        ¥ ( value -- )  Fills all elements with value.
  730.     idxbase  limit sfloats  bounds
  731.     ?DO  fdup  i f!  4 +LOOP  fdrop  ;m
  732.  
  733. :m PRINT:    ¥ Prints all elements
  734.     limit: self 0 ?DO  i dup  4 .r  space  at: self  e. cr
  735.     LOOP  ;m
  736.  
  737. :m CLASSINIT:
  738.     undef
  739.     limit: self  FOR  fdup  i to: self  NEXT  fdrop  ;m
  740.  
  741. ;class
  742.  
  743.  
  744. ¥ Fmatrix implements a 2-dimensional FP matrix.
  745. ¥  Usage:   50 100  dimension   Fmatrix MM
  746.  
  747.     0    value    ROWDIM
  748.     0    value    COLDIM
  749.  
  750. : DIMENSION
  751.     -> colDim   -> rowDim
  752.     colDim  rowDim  *  ;
  753.  
  754.  
  755. :class    FMATRIX    super{ Farray }
  756.  
  757.     var    #rows
  758.     var    #cols
  759.     var    rowLength
  760.     var    colLength
  761.  
  762. :m #ROWS:    get: #rows   ;m
  763. :m #COLS:    get: #cols   ;m
  764.  
  765. :m ^ELEM:  { row col ¥ temp -- addr }
  766.     row  get: #cols  *  col +  ^elem  ;m
  767.  
  768. :m AT:    ^elem: self  f@  ;m
  769. :m TO:    ^elem: self  f!  ;m
  770.  
  771. :m ROW:  { row -- limit start stride }    ¥ Sets up for a DO over the row.
  772.     row  get: rowLength *  idxBase +        ¥ addr
  773.     get: rowLength  bounds                    ¥ ( limit addr )
  774.     1float   ;m                                ¥ stride
  775.  
  776. :m COL:  { col -- limit start stride }    ¥ Sets up for a DO over the column.
  777.     col floats  idxBase +                    ¥ addr
  778.     get: colLength  bounds                    ¥ ( limit addr )
  779.     get: rowLength  ;m                        ¥ stride
  780.  
  781. :m PUTROW:    ¥ ( 1st ... last  row -- )
  782.     row: self   drop  ( we know stride is 1float )
  783.     swap  1float -
  784.     DO  i f!   1float negate  +LOOP  ;m
  785.  
  786. :m PRINTROW:  { row -- }
  787.     row  row: self    drop
  788.     ?DO   i f@ 10 e.r  1float +LOOP    ;m
  789.  
  790. :m PRINTCOL:  { col ¥ stride -- }
  791.     col  col: self  -> stride
  792.     ?DO  i f@  10 e.r  cr  stride +LOOP  ;m
  793.  
  794. :m PRINT:
  795.     get: #rows  0
  796.     ?DO  i printRow: self  cr  LOOP  ;m
  797.  
  798. :m CLASSINIT:
  799.     rowDim  put: #rows   colDim  put: #cols
  800.     get: #cols floats  put: rowLength
  801.     get: #rows  get: rowLength *  put: colLength
  802.     classinit: super   ;m
  803.  
  804. ;class
  805.  
  806.  
  807. endload
  808.  
  809.  
  810. 3 5 dimension  fmatrix  matA
  811. 5 3 dimension  fmatrix  matB
  812. 3 3 dimension  fmatrix  matC
  813.  
  814. -2.0 -7.0 -4.0 8.0 1.0      0  putRow: matA
  815. -3.0 0.0 -1.0 0.0 -9.0        1  putRow: matA
  816. 0.0 -1.0 -8.0 -9.0 -3.0        2  putRow: matA
  817.  
  818. -1.0 3.0 -2.0                0  putRow: matB
  819. 4.0 -1.0 1.0                1  putRow: matB
  820. 9.0 -4.0 -8.0                2  putRow: matB
  821. -5.0 0.0 -3.0                3  putRow: matB
  822. 6.0 3.0 -6.0                4  putRow: matB
  823.